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We have studied transition metal clusters from a quantum information theory perspective us- 
ing the density-matrix renormalization group (DMRG) method. We demonstrate the competition 
between entanglement and interaction localization. We also discuss the application of the configura- 
tion interaction based dynamically extended active space procedure which significantly reduces the 
effective system size and accelerates the speed of convergence for complicated molecular electronic 
structures to a great extent. Our results indicate the importance of taking entanglement among 
molecular orbitals into account in order to devise an optimal orbital ordering and carry out efficient 
calculations on transition metal clusters. We propose a recipe to perform DMRG calculations in 
a black-box fashion and we point out the connections of our work to other tensor network state 
approaches. 



I. INTRODUCTION 

The most important property of multi-component 
quantum systems is entanglement which corresponds to 
quantum correlations between particles or a collection 
of particles forming a larger subsystem The degree 
of entanglement is decisive for the behavior of all multi- 
component quantum systems and for the numerical algo- 
rithms developed to simulate such systems 

The concepts of quantum information theory (QIT) 
have been introduced to the density matrix renormaliza- 
tion group (DMRG) method 0, y] and created a fresh 
impetus to the development of new methods that focus 
on entanglement optimization (l2l.[l3|. It was shown that 
the quantum information entropy is a direct measure of 
the behavior of quantum systems. Thus, it is mandatory 
to develop methods and techniques to obtain as much 
information as possible from entropy analyses [l3l - [l8j . 

Since the first application of the DMRG algorithm 
to quantum chemical systems using the full electronic 
Hamiltonian (l9| . the method has gone through major 
algorithmic developments by various groups [20|-[2J]. In 
2002, two groups independently provided the first large 
scale calculations on diatomic molecules [25|, [26[ . Calcu- 
lations on the F2 molecule by keeping more than 2000 
block states (i.e., states defined on the active subsystem 
in DMRG) denoted as M in the following, and on the wa- 
ter molecule using 6000 block states [27| confirmed that 
the DMRG method is capable of reaching the limit of 
the full-CI method in cases where most standard quan- 
tum chemical approaches fail and cannot be applied due 
to the requirement of very large active spaces. Effi- 
cient calculations of excited states and the relation be- 
tween the DMRG wave function and the standard CI- 
expansion have already been discussed in Ref. [28j |. The 
authors have also shown that the DMRG wave function 
is also suited to study problems when the characteris- 
tics of the wave function changes drastically, for exam- 



ple, when the bond length between diatomic molecules is 
changed. Therefore, it became evident that very accurate 
results can be obtained by increasing the computational 
resources related to the DMRG block states. On the 
other hand, the more important questions, namely what 
is the minimum computational effort to obtain results 
with a given accuracy is still unsolved and part of active 
research [H, GJ, [H, Q . 

A key-ingredient of the DMRG method related to its 
multi-component subsystem wave function representa- 
tion is entanglement. A controlled manipulation of it 
is expected to extend capabilities of DMRG to treat 
larger systems in a more efficient way. The develop- 
ment of entanglement reduction algorithms (ERA) is 
thus very appealing. Matrix product states (MPS) [3(3] 
that are inherently produced by the DMRG algorithm 
can be used to localize the entanglement by reordering 
sites lEIJF 
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25. l26l. 311 or by optimization of the ba- 
j32-37|7^ Considering the entanglement 
between arbitrarily chosen pairs of sites [l^, HU , a net- 
work topology can be determined for tensor-network- 
state (TNS) algorithms [Tl|. Further gain in speed of 
convergence can be achieve d by an optimization of the 
initialization procedure (T3 , l38j in which highly entan- 
gled states are taken into account from the very begin- 
ning. One such algorithm is the dynamically extended 
active space (DEAS) method and its extensio n by 
including Cl-expansion techniques (CI-DEAS) [H, flof . 
In this article, we will show that the proposed entropy- 
based optimization scheme including the (CI-)DEAS and 
the dynamic block state selection (DBSS) procedures are 
vital ingredient to obtain highly accurate results with 
tremendous savings in computation resources and time. 
In addition, it is a smart way of reducing the Hilbert 
space which allows us to study large active spaces with 
a modest number of renormalized states in contrast to 
inefficient and expensive brute-force DMRG calculations 
at a fixed number of renormalized states. DMRG cal- 
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culations on complex chemical systems using quantum 
information theory have not yet been carried out, and 
are the subject of this work. 

The very challenging binuclear oxo-bridged copper 
clusters identified in Ref. as a very difficult case for 
complete-active-space-like calculations pose an ideal test 
for the DMRG algorithm. In a pioneering study, some of 
us investigated the DMRG algorithm for the prediction 
of relative energies of transition metal clusters of different 
molecular structures, the bis[u-oxo) and fi — rj 2 : rj 2 per- 
oxo isomers of [Cu2 02] 2+ [42|]. Based on these results, 
we drew the conclusion that the DMRG algorithm is 
suited to study transition metal complexes and clusters. 
Kurashige and Yanai have picked up that same problem 
as well and performed massively parallel DMRG calcula- 
tions employing 2400 renormalized basis states but on a 
larger active space with a smaller one-particle basis [43[ . 
In their study, a new aspect was the inclusion of a per- 



turbative correction introduced originally by White 44] 
for the one-site DMRG algorithm to prevent convergence 
to local minima during the optimization process. In the 
present study of [Cu2 02] 2+ , we show that no convergence 
acceleration technique such as white noise (25j or per- 
turbative correction [43| is needed if an entropy-based 
optimization scheme is employed. 

The purpose of the paper is to show how quantum 
information entropies can be efficiently exploited to per- 
form accurate calculations on complicated chemical sys- 
tems such as transition metal clusters. Our results in- 
dicate that - within chemical accuracy - much larger 
molecular systems can be studied than before using very 
little computational resources. In Sec. Ill, we briefly 
discuss the various technical aspects of our calculations 
while the entropy analysis of the two isomers are pre- 
sented in Sec. IV. In Sec. V, extensions to tensor net- 
work state methods are discussed and our conclusions 
are given. 



II. THE CHALLENGE OF BINUCLEAR 
COPPER CLUSTERS 

The reliable first-principles description of transition 
metal complexes and clusters remains an important task 
for theoretical molecular physics and quantum chemistry 
- especially because of the role of such molecules in 
catalysis and bioinorganic chemistry [45l - |48j . 

In 2006, Cramer et at investigated several theoret- 
ical models on the [Cu2 02] 2+ torture track and found 
incisive discrepancies between them [4l|. The CASSCF 
calculations even yield a qualitatively wrong interpreta- 
tion of the energy difference between the two isomers. 
The reason for this striking failure lies in the inability 
of the CASSCF method to include all relevant molecular 
orbitals into the one-particle active space that would be 
necessary to obtain a qualitatively correct description of 
the electronic structure. Already the qualitative picture 
of an extended Hiickel calculation indicated the require- 



ment of more than doubling the active space for binuclear 
transition metal clusters compared to the mononuclear 
analog [42j. An additional support for the fact that a 
very large active space is needed for a correct descrip- 
tion of the electronic structure of transition metals is the 
finding of Pierloot et al. to include a second d shell to ob- 
tain accurate results (49j |. As we showed in our previous 
work [iij], the DMRG algorithm is an ideal candidate to 
tackle the description of challenging electronic structures 
that require large active spaces occurring in transition 
metal chemistry. 

This work reinvestigates the dicopper clusters with 
an improved methodology transforming the DMRG ap- 
proach eventually into a black-box method with signif- 
icantly reduced computational requirements. For our 
DMRG calculations, the same active space, one-particle 
basis set and effective core potential as in Ref. [42j was 
used. All quantities calculated in this paper will be given 
in Hartree atomic units of which the energy unit is one 
Hartree. 



III. NUMERICAL PROCEDURE 

In this section, we outline the procedures and meth- 
ods needed for the entropy-based optimization scheme to 
efficiently carry out DMRG calculations. Examples and 
figures will be presented for the /i — i] 2 : rj 2 peroxo iso- 
mer of [Cu202] 2+ while for the bis(^t-oxo) isomer only 
the final results are given. 



A. Molecular Hamiltonian in second quantization 

In quantum chemical (QC-)DMRG applications, the 
electron-electron repulsion is taken into account by an 
iterative procedure that minimizes the Rayleigh quotient 
corresponding to the electronic Hamiltonian given by 
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and thus determines the full-CI wave function. In Eq. ([1}, 
Ci a and c\ a are the usual electron annihilation and cre- 
ation operators, Ty denotes the one-electron integral 
comprising the kinetic energy of the electrons and the 
external electric field of the nuclei. Viju stands for the 
two-electron integrals and contains the electron-electron 
repulsion operator, defined as 



, , 3 ,3 <f>i($i)<tij(%2)Mx2)Mxi) 
Vijki = I <fxid A x 2 - 



\xi - x 2 \ 



(2) 



We obtain the Hartree-Fock orbitals in a given one- 
particle basis of atomic Gaussian functions and trans- 
form the one-electron and two-electron integrals in the 
atomic basis set to the Hartree-Fock molecular orbital 
basis using the MoLPRO program package [SOj . which 
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FIG. 1: (Color online) Pictorial representation of the magni- 
tude of components Ty of the one-electron operator for the 
fi — rj 1 : rj 2 peroxo isomer of [Cu2C>2] 2+ . For better visibility 
the one-dimensional orbital chain is plotted on a circle with 
modulated radius in a clockwise direction. Orbitals belonging 
to different irreducible representations are shown by different 
symbols. Numbers next to the symbols label molecular or- 
bitals. 

we also employ to obtain reference complete-active-space 
self-consistent field (CASSCF) and complete-active-space 
configuration-interaction (CASCI) energies. 

In the QC-DMRG algorithm, a one-dimensional chain 
is built up from molecular orbitals obtained from a suit- 
able mean-field or multi-configuration self-consistent field 
(MCSCF) calculation. As will be discussed later, the one- 
orbital entro py function [l3| and the two-orbital mutual 
information [1 51 ] provide a good starting configuration. 
The irreducible representations of the orbitals can also 
be used in the DMRG procedure to carry out calcula- 
tions for a given point-group symmetry |13| . This will be 
used in the present work, and as an example, the compo- 
nents Tij of the one-electron operator are shown in Fig.Q] 
for the n — rj 2 : rj 2 peroxo isomer of [Cu202] 2+ . 



B. Two-site variant of the DMRG method 

For our DMRG calculations, we employ the two-site 
variant of the DMRG method introduced by White ||. 
In the superblock configuration, two sites are between 
two blocks (sets) of orbitals. To fix the notation for the 
rest of the paper, we label the block states of the left 
block containing I orbitals, the two intermediate sites and 
the right block with r orbitals by Mi, qi+i, qi+2, M r , re- 
spectively. The total number of orbitals N is I + 2 + r. 
When a larger block is formed during the rcnormaliza- 
tion steps by adding a single site to the left or to the 
right block the new block states are denoted by Ml and 
Mr where L = I + 1 and R = r + 1. The one-site vari- 



ant of the DMRG [3 related to the MPS approach was 
first applied in quantum chemistry in Ref. [261 ] for fixed 
number of block states but convergence properties were 
found considerably worse than for the two-site DMRG al- 
gorithm. This can be improved by the introduction of a 
local density operator [44j which was used by Kurashige 
and Yanai for the calculation of the transition metal clus- 
ter compound with 2400 block states [43j. Although the 
size of superblock Hilbert space is smaller in this case, we 
use the two-site variant since this configuration allows 
one to control the number of block states dynamically 
and convergence to local minima is less likely f26j . 

C. Dynamic Block State Selection 

A fundamental concept related to the inseparability 
and non-locality of quantum mechanics is entanglement. 
Since the QC-DMRG algorithm approximates a compos- 
ite system with long-range interactions, the results of 
quantum information theory can be used to understand 
the criteria of its convergence. 

The two-site variant of the DMRG method has origi- 
nally been employed with a fixed number of block states 
while the degree of entanglement between the DMRG 
blocks for a given superblock configuration is related 
to the Schmidt number rg c hmidt- For a pure target 
state |*ts> e A = A L (g> A R , with dimA L = M L , 
dimA/j = Mr, the Schmidt decomposition states that 

Schmidt <min(M L ,M H ) 

I*ts}= Yl Wi|ei>®l/*), (3) 

i=i 

where |e») <8> \fi) form a bi-orthogonal basis (e-i\ej) — 
(filfj) — fiiji an( l < 0;^ < 1 with the condition 
J2iUi 2 = 1. If r Sc hmidt > 1 then according to Ref. [51] 
|\&ts) is inseparable and the two blocks are entangled. 

Since possible measures of entanglement for fcrmionic 
systems are the von Neumann and Renyi entropies, it 
is more efficient to control the truncation error (2(| or 
the quan tum information loss x a t each renormalization 
step [Tj] . In the DMRG procedure, during the renormal- 
ization step the block Bl is formed of the subblock Bi 
and the I + 1 th site. Denoting by sl(1) the entropy of the 
left subblock of length I and by s; + i the entropy of the 
I + 1 th site, the change of the block entropy by forming 
a larger block, Bl(1 + 1), is given as 

s L (l) + s(l) i+1 + I L (l) = s L (l + 1) , (4) 

where the von Neumann entropy of a block with I orbitals 
is given as 

slQ) = -^w„lnw a , (5) 

a 

where u> a stands for the eigenvalues of the reduced den- 
sity matrix of the block Bi . The so-called mutual infor- 
mation < quantifies the correlation between the 
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subsystem and the site and it is zero if and only if the 
two blocks are uncorrelated. 

In order to control the quantum information loss, Ml 
(or Mr) is increased systematically at each renormaliza- 
tion step until the following condition holds 



s L (l + l)- S l*™ c (l + l)<x, 



(0) 



where \ is an a priori defined error margin. For + 1), 
i.e., before the truncation, a — 1 . . .Mi x qi+\ while for 



sf unc (l + 1) according to Eq. © M^ unc < Mi x q l+1 
is used. This approach guarantees that the number of 
block states are adjusted according to the entanglement 
between the DMRG blocks and the a priori defined ac- 
curacy can be reached 
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In addition, an entropy sum 
rule can be used as an alternative test of convergence [3] . 
Therefore, we set the minimum number of block states 
to M m i„ and x- By setting M m i D ~ q 3 or q 4 , convergence 
to local minima can be avoided. In our implementation, 
we have q = 4 and the basis states correspond to the 
1 0) , | 4-), t) an d I 4-T) states. The maximum number 
of block states selected dynamically during the course 
of iterations will be denoted by M max . This approach, 
however, does not work for the one-site variant of the 
DMRG algorithm since the Schmidt number of a one- 
site superblock configuration Ml = Mi x qi +i cannot be 
larger than M r . This prevents Mi to increase above M r 
according to Eq. ©. 



D. Entanglement and interaction localization 

The von Neumann and Renyi entropies have been stud- 
ied for quantum chemical systems as well and it was 
shown that orbitals lying closer to and further away from 
the Fermi surface possess larger and smaller orbital en- 
tropy, respectively !13]. The orbital entropy is related 
to the mixture of a local state and it is expressed by 
the eigenvalues of the reduced density matrix for a given 
orbital, namely 



s(l)j = ~2^UW 



In i 



(7) 



where i = 1 . . . N labels the orbital index while uj a> i 
stands for a = 1 . . . qt the eigenvalues of the reduced 
density matrix of orbital i. In Fig. [21 the single-orbital 
entropy is shown for the /i — rj 2 : rj 2 peroxo isomer of 
[Cu2C>2] 2+ calculated by setting the quantum informa- 
tion loss x = 10 -5 . Orbitals with a large entropy sig- 
nificantly contribute to the correlation energy whereas 
other slightly entangled orbitals do not. Since the or- 
bitals possess different single-orbital entropies, the or- 
dering of orbitals along the one-dimensional chain of or- 
bitals in the DMRG algorithm has a major impact on 
the block entropy s(l), so that the block entropy profile 
can be changed based on the ordering of the orbitals [l3[ . 
The block entropy also determines the required compu- 
tational resources to reach a given accuracy [l3|, [l4| • If 
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FIG. 2: (Color online) Single-orbital entropy function ob- 
tained for the (i — rj 2 : rj 2 peroxo isomer of [Cu2C>2] 2+ by 
setting the quantum information loss \ — 10 -5 . The symbols 
indicate corresponding point-group symmetries of the ener- 
getically ordered orbitals. 



an optimized ordering is used, DMRG results can be 
obtained using a considerably smaller number of block 
states for a given error bound. 

An optimal orbital alignment by means of speed of con- 
vergence can be obtained by reordering the orbitals, so 
that the DMRG blocks are entangled only for a few iter- 
ations. This can be achieved by placing highly entangled 
orbitals at one of the ends or close to the center of the 
chain. 

Unfortunately, this is not true in general, since the 
independent interaction terms, like Ty and the local di- 
rect, pair-hopping, and spin-flip terms of Vijki act as in- 
dependent quantum channels; they all generate different 
amounts of entanglement [1 61 ] - Hence, localizing entan- 
glement generated by one channel, e.g., the one-electron 
term, might lead to delocalized entanglement in another 
channel. As an example, the one-electron operators are 
analyzed for the fi — rj 2 : rj 2 peroxo isomer of [Cu202] 2+ 
shown in Fig. [TJ In this case they do not couple orbitals 
among different irreducible representations of the D2h 
point group while the two-electron Coulomb repulsion 
operator does. A proper cost function to take care of in- 
teraction and entanglement localization can be expressed 
based on the two-orbital mutual information, 



hj = s(2) itj - s(l)i - 



j • 



(8) 



where s(2)ij is the two-orbital entropy between a pair of 
sites which was introduced to the QC-DMRG by Rissler 
et al. [15]. If the electron-electron interactions are ne- 
glected, the two-orbital mutual information has a similar 
structure as shown in Fig. [TJ where only the orbitals of 
the same irreducible representations are entangled. The 



5 



O 



(a) 



5 10 15 20 25 30 35 40 45 

Orbital index 





(b) 






v/*nA/vJ 






i fr-£rw '-w 1 — 








FIG. 3: (Color online) (a) The single-orbital entropy, (b) 
block-entropy, and (c) the mutual information obtained for 
the (i — rj 1 : rf peroxo isomer of [Cu2C>2] 2+ after the seventh 
DMRG sweep while neglecting the electron-electron interac- 
tions. 



resulting single-orbital entropy, block-entropy, and the 
mutual information are shown in Fig. |3] obtained by the 
DMRG method after the seventh sweep. The block en- 
tropy oscillates and is exactly zero when all orbitals of 
the same irreducible representations belong to the left 
or to the right DMRG block. Therefore, the total wave 
function can be expressed as a product state of the wave 
functions of the subblocks of irreducible representations. 
When the two-electron integrals are also considered, or- 
bitals among different irreducible representations are also 
entangled as shown by the components of the two-orbital 
mutual information in Fig. U for the \i — r/ 2 : rj 2 peroxo 
isomer of [Cu202] 2+ . Since the one-electron integrals are 
usually an order of magnitude larger in chemical sys- 
tems, an optimal ordering can be found by reordering 
the orbitals within the same irreducible representation 
and reorderin g th e blocks of different irreducible repre- 
sentations way, the block entropy can 
be reduced significantly but the one-electron terms hav- 
ing the strongest interaction remain as local as possible. 
To draw the analogy to the Hubbard model, most chem- 
ical systems would therefore correspond to the so-called 
small U limit. 

The reordering concept can be put in a more rigorous 
form by minimizing the entanglement distance, which can 
be expressed as a cost function, 



Jdist — 



(9) 
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where the entanglement between pair of sites is weighted 
by the distance in the chain between the orbitals. In 



FIG. 4: (Color online) Components of the two-orbital mutual 
information which are larger than 10 -4 for the fj, — r/ 2 : r/ 2 
peroxo isomer of [Cu202] 2+ obtained with x = 10 -5 and for 
the energetical ordering. The role of symbols and numbers is 
the same as in Fig. [1] 



Ref. [Tf|, the effect of the parameter 77 = — 2 has been an- 
alyzed using simulated annealing. Additional parameters 
to weight the off-diagonal elements of I have been studied 
as well. In our approach, we use both rj = 1 and 2. The 
latter choice also has the advantage that it can be related 
to the spectral algorithms of seriation problems [52| . The 
main aim is to sequence a set of N objects, i.e., to bijec- 
tively map the elements to the integers 1, . . . ,N based 
on a real valued correlation function f(i,j) = f(J,i) 
which reflects the desire for items i and j to be near 
each other in the sequence. The two-orbital mutual in- 
formation is such a function which can also be seen as 
a weighted graph. In general, the problem of finding all 
ways to sequence the elements, so that the correlations 
are consistent, becomes NP-complete [531 ] in the presence 
of inconsistencies. In such a case, there may be no con- 
sistent solution and one needs to find the best approxi- 
mation. If a consistent ordering is possible, the problem 
is well posed. Most of the combinatorial algorithms for 
well-posed problems break down when the data is incon- 
sistent, limiting their value for many problems. In our 
approach, the minimization is performed iteratively with 
the constraint that orbitals of the same irreducible repre- 
sentations are kept together, thus reordering of orbitals 
is allowed within irreducible representations and reorder- 
ing the blocks of orbitals of the irreducible representation 
is also allowed. As a result, a highly optimized ordering 
can be obtained as is shown in Fig. [5] The value for the 
energetical ordering idist = 821.4 is reduced to 134.1 us- 
ing rj — 2. A smaller value of 63 could be reached by 
excluding the constraint discussed above. However, the 
DMRG calculations perform considerably worse for that 
ordering. 

A further justification of our cost function can be dis- 
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entropy profile for a given ordering of molecular orbitals 
can be determined as a function of sweeps [HI, |39| . Once 
the wave function has converged by means of the entropy 
sum rule (l4| . the single-orbital entropy profile for the 
given target state can be obtained. This well known pro- 
cedure of collecting data points from subsequent renor- 
malization steps of a full sweep was also used by Ghosh 
et al. to efficiently calculate four-point correlation func- 
tions [37]. The one-orbital entropy, on the other hand, 
can also be expressed in terms of the occupation-number 
representation [l5j]. Therefore, the calculation of s(l)j is 
also possible once the required operators are determined 
for the given superblock configuration. A direct com- 
parison of data points obtained by the two approaches 
provides a reliable error estimate. 



FIG. 5: (Color online) Similar to Fig. [4] but for the reordered 
orbitals. 

cussed in terms of spectral algorithms. Since the mini- 
mization of the cost function g(n) — J2ij fihj)( n i~ 7r j) 2 
is hard in terms of permutations 7T, it can be approx- 
imated by a cost function like /^ist of continuous vari- 
ables Xi that maintains its structure. From spectral 
graph theory it follows that the so called Fiedler vector 
x = (xi, . . .xn) is the solution that minimizes F(x) — 
x' Lx — y~]{j Ii,j{xi — Xj) 2 subject to the following con- 
straints that . Xi = and x i = 1> where the graph 
Laplacian is L = D — I with — J^j h,j- The second 
eigenvector of the Laplacian is the Fiedler vector [EH, l55j 
which defines a (I-dimensional) embedding of the graph 
on a line that tries to respect the highest entries of Iij 
and the edge length of the graph. Ordering the entries 
of the Fiedler vector by non-decreasing or non-increasing 
way provides us a possible ordering. A naive applica- 
tion of optimization methods based on the Fiedler vector 
yielded a worse ordering than the one shown in Fig. [5] 
A more detailed analysis of energetical ordering based on 
the Fiedler vector will be part of our subsequent work. 
Numerical results obtained by the optimized ordering will 
be further discussed in Sec. IV. 

Note that a new electronic wave function ansatz in 
terms of the complete-graph tensor network (CGTN) 
parametrization [10] might contain the same information 
as the weighted graph of the cost function in Eq. ([9]). 
This efficient parametrization might also be used to de- 
vise an optimized orbital ordering. 



E. Efficient calculation of the single-orbital entropy 

The single-orbital entropy can be calculated for each 
renormalization step of a full sweep, thus s(l)i can be 
obtained for i = 1...N. Therefore, the single-orbital 



F. Efficient calculation of the two-orbital entropy 
function 

The two-orbital entropy function can be expressed in 
terms of the occupation- number representation [l5j . Un- 
fortunately, it requires the calculation of 23 independent 
two-point correlation functions. Since all two-point cor- 
relation functions have to be renormalized and stored 
in memory or on disk due to the truncation of the su- 
perblock Hilbert space, the efficient calculation of these 
operators is crucial for feasible calculations with respect 
to wall time and computational resources. 

In a standard real space DMRG procedure, the cor- 
relation functions are usually calculated for the sym- 
metric superblock configuration, i.e., when the size of 
the left and right blocks are equal. This configuration 
provides the most accurate result for a fixed number of 
block states since the block entropy reaches its maximum 
value, so that the highest level of entanglement can be 
reached [3]. In fact, the largest error in a measurable 
quantity is related to the largest truncation error within 
a full DMRG sweep [56|. In contrast to this, using the 
DBSS method, the error can be kept below an a priori de- 
fined threshold, and hence an accurate calculation of the 
correlation functions is possible for the non-symmetric 
superblock configurations as well. In addition, if two non- 
interacting orbitals are attached to both ends of the chain 
the reduced block density matrix at the turning points 
of a sweep, where all orbitals belong to the left or the 
right block, has only one non-zero eigenvalue according 
to the Schmidt decomposition. Renormalized operators 
required for the two-orbital entropy function reduce to 
single numbers (M; = M r — 1) [Til, Hs|j consequently. 
Even without this trick, an efficient calculation of the 
correlation function at the turning points is possible be- 
cause the environment block contains q basis states and 
the system block only M min states when the DBSS proce- 
dure is used. This is crucial in the quantum information 
analysis since the construction of the mutual entropy is 
very costly. In our approach, all the required operators 
are generated only for the superblock configuration at the 
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FIG. 6: (Color online) Number of entanglement bonds emerg- 
ing from the orbitals based on Fig. [4] The meaning of the 
symbols is the same as in Fig. [2] 



turning points for which the correlation functions are cal- 
culated. This is achieved by subsequent renormalization 
of the operators based on the a priori determined trans- 
formation matrices 301. 



G. Configuration Interaction based Dynamically 
Extended Active Space (CI-DEAS) procedure 

The non-local version of the DMRG method is very 
sensitive to the initialization procedure. If a poorly ap- 
proximated starting configuration is used, the conver- 
gence can be very slow and the DMRG can even be 
trapped in local minima [38j |. This can, however, be 
avoided by including highly entangled orbitals from the 
very beginning and expanding the active space dynami- 
cally [13|. This approach has also been extended to in- 
clude protocols based on the Configuration Interaction 
(CI) procedure [39|, H(| and applied to systems with sizes 
up to 60 orbitals. Here, we briefly summarize the main 
aspects of the method as required for the discussion in 
this work. 

Taking a look at Fig. @J some orbitals are highly entan- 
gled with several other orbitals while others are entan- 
gled with a few orbitals only. The number of entangle- 
ment bonds emerging from the orbitals based on Fig. @] 
is shown in Fig. [5] 

In order to guarantee fast convergence the highly 
entangled orbitals — those with several entanglement 
bonds — should be included from the very beginning 
of the calculations. If the bonds are also weighted with 
their strength, by summing / column-wise, one gets back 
the diagonal entries of the graph Laplacian which corre- 
sponds to the single-orbital entropy. In the DEAS proce- 



dure, we introduce a complete-active-space (CAS) vector 
that selects the highly entangled orbitals. The CAS vec- 
tor (CASV) is formed by ordering orbitals with decreas- 
ing entropy values. 

Since the DMRG algorithm itself is a basis-state trans- 
formation protocol, that transforms single-particle ba- 
sis states to mult i- particle basis states, the environment 
block in the DEAS procedure is formed for each DMRG 
iteration step from the one-particle basis states of those 
orbitals which posses the largest single-orbital entropies. 
In the first iteration step the left block (system block) 
contains one orbital represented by q states while the 
right block (environment block) r = N — I — 2 orbitals. 
Since the exact representation of the right block would 
require q N ~ l ~ 2 states which is too large for large sys- 
tems, only a subset of orbitals is included to form the 
active space. Therefore, at each DMRG iteration step of 
the warm up procedure, i.e., for iterations 1 ... N — 3 the 
M r states are formed from those components of the CAS 
vector which belong to the environment block and pos- 
sess largest entropies. The starting value of M r (M sta rt) 
is set prior of the calculation but during the iteration pro- 
cedure M r is adjusted as M r = max(M/, M star t) in order 
to construct at least as many environment states as the 
left block has in order to satisfy the constraints set by 
the Schmidt decomposition. Identifying orbitals of the 
right block as doubly filled, empty or active orbitals the 
effective size of the environment block can be reduced 
significantly. The empty orbitals can be neglected, while 
a partial summation over the doubly filled orbitals results 
in a constant shift in the energy. If DMRG auxiliary op- 
erators [57j j are formed by partial summations on the left 
block, the effective system size of the problem is deter- 
mined by the active orbitals only [l3|, [33, S3 ■ Therefore, 
in the warm up procedure the effective size of the sys- 
tem is reduced to 5 to 7 orbitals which allows one to use 
larger Af sta rt without a significant increase in computa- 
tional time. 

In order to construct even better environment states, 
we also utilize CI expansion procedures. In standard CI 
techniques, the trial wave function is written as a linear 
combination of determinants with expansion coefficients 
determined by requiring that the energy should be mini- 
mized [581 ] . The molecular orbitals used for building the 
excited Slater determinants are taken from a CASSCF 
calculation and kept fix. Therefore, in the CI method 
the number of determinants included in the treatment 
is increased systematically in order to achieve a better 
accuracy. Determinants which are singly, doubly, triply, 
quadruply, etc. excited relative to the Hartree-Fock (HF) 
configuration are indicated by the subscripts S, D, T, Q. 
The exact wave function in a given one-particle basis, the 
full- CI wave function, is then given as 

*ci = oo*hf + ^2 as ® s + aD ^ D + ^ t ^ t + ■■■ 

S D T 

(10) 

Since the segment of the HF-orbitals belonging to the 
environment block is known, the restricted subspace of 
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the environment block can be formed for a given Cl-level 
in the DEAS procedure. Therefore, the right block con- 
tains states for a given Cl-level while the total wave func- 
tion can contain higher excitations as well due to the cor- 
relation between the two blocks. The environment block 
states are constructed at each iteration step, so that the 
environment block is always optimized for the renormal- 
ized system (left) block. This procedure guarantees that 
several highly entangled orbitals are correlated from the 
very beginning and both static and dynamic correlations 
are taken into account. The reduced density matrix is 
well defined, thus block states can be selected efficiently 
based on the entropy considerations and convergence to 
local minima can be avoided. Since a significant part 
of the correlation energy can be obtained in this way, 
usually at the end of the initialization procedure, i.e., af- 
ter one-half sweep, chemical accuracy is often reached. 
The initial CAS vector can be determined based on the 
chemical character of the molecule or in a self-consistent 
fashion. In the latter approach, the CAS vector is set first 
as CASV = [N,N—1,...,1] to include long range inter- 
actions between the left block and the rightmost sites of 
the chain and a calculation using small values of M m j n 
and Mstart is performed. After a full sweep, the entropy 
functions are calculated and the ordering as well as the 
CAS vector are determined. Even though the DMRG 
wave function has not yet converged, most relevant in- 
formation of the system can already be extracted. In 
addition, the DMRG results can be systematically im- 
proved using the optimized ordering and the CAS vector 
as a starting point for new DMRG calculations. The CI- 
DEAS method allows a simple and fast calculation of all 
physical quantities at the end of the first sweep. 

In the following section, the DBSS and CI-DEAS pro- 
cedures are applied to chemical systems and the results 
are discussed in more detail. 



IV. RESULTS 

In this section, we present results for the two isomers of 
[Cu 2 2 ] 2+ . We discuss our DBSS/CI-DEAS procedures 
by analyzing the peroxo isomer while for the bis(/^-oxo) 
isomer only the final results will be given. 

A. Electronic ground state of the two [Cu2C>2] 2+ 
cores 

Here, we summarize the main steps of our entropy- 
based optimization procedures. First, we perform a short 
DMRG calculation using a limited set of block states in 
order to obtain the most important characteristics of the 
entropy functions and to determine the low-lying energy 
spectrum. In order to use DMRG in an automated man- 
ner as a black-box method, we do not use any specific 
ordering or CAS vector related to the chemical charac- 
teristics of the problem. 
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FIG. 7: (Color online) (a) Convergence of the ground state 
energy in Hartree for various parameters for the fx — r/ 2 : r/ 2 
peroxo isomer. The dashed line gives the ground state energy 
obtained in Ref. 42] using M = 800 block states, (b) Number 
of block states as a function of iteration steps for % — 10 -4 
corresponding to Fig. [9] 

In order to determine the multiplicity of the converged 
target state, we also calculate the expectation value of 
the S 2 operator given as 

s 2 = + s;s! (ii) 

ij ij i 

The expectation value is equal to S(S + 1) in Hartree 
atomic units, i.e., zero for a singlet state and two for a 
triplet state. 

We performed a calculation with fixed M = 64 states, 
using the energetical ordering, and the CAS vector was 
simply set to CASV= [N, N-l, . . . , 1]. This latter choice 
guarantees that during the initialization procedure long 
range interaction between the system block (left-block) 
and the rightmost sites of the chain is taken into account. 
In this calculation, we have restricted the CI-DEAS pro- 
cedure to include determinants with at most triple exci- 
tations from the Hartree-Fock state. The convergence of 
the ground state energy in the A g irreducible representa- 
tion of the point group D2h is shown in Fig. EJa) as the 
square symbols. Although the convergence is rather slow, 
the results obtained in our previous work with M = 800 
block states [42| - shown by the dashed line - can be 
reached after eight sweeps. The convergence to the sin- 
glet state has been confirmed by the expectation value of 
the total-spin operator (S 2 ) — 10~ 3 . 

The block entropy profile converges very slowly as 
shown in Fig. [5Jb). In fact, the entropy profile cor- 
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FIG. 8: (Color online) (a) Single-orbital entropy, (b) block en- 
tropy and (c) mutual information profile for the fi — rf : rf 
peroxo isomer of [Cu2 02] 2+ using energetical ordering, with- 
out optimized CAS vector for the singlet ground state of A g 
symmetry with fixed M — 64 block states. The convergence 
of the ground state energy in Hartree is shown in Fig. [TJa) by 
the square symbol. 



FIG. 9: (Color online) (a) Single-orbital entropy, (b) block- 
entropy, and (c) mutual information obtained for the fj, — -q 2 : 
rj 2 peroxo isomer of [Cu 2 2 ] 2+ as a function of DMRG sweeps 
with optimized ordering, CAS vector, and by setting the 
quantum information loss % — 10 -4 . The CI-DEAS initializa- 
tion procedure corresponds to symbols with red circle. The 
optimized ordering and CAS vector is given in Table [U 



responding to the CI-DEAS procedure represented by 
the red circle symbol is almost zero for most of the su- 
perblock configurations which indicates the lack of entan- 
glement between the DMRG blocks. As a consequence, 
the DMRG algorithm selected inappropriate states for 
the description of the environment blocks. 

A detailed analysis of the entropy functions provides 
even better convergence properties. The obtained single- 
orbital and block entropies as a function of DMRG 
sweeps are shown in Fig. [8] By comparing Figs.[H[a) and 
[3J we see that orbitals with large entropies can already 
be identified and a reasonable CAS vector can thus be 
constructed. In Fig. Wp), one also recognizes that once 
the third orbital is pushed in the left block, i.e., for I > 3 
the block entropy increases significantly. This is because 
the third orbital is highly entangled with the 14 th and 
25 th orbitals as can be seen in Fig. 0J The largest values 
of the block entropy is reached for 13 < I < 35 when 
the highly entangled orbitals 13 and 35 belong to the 
two different blocks. In order to localize entanglement, 
ordering optimization based on the two-orbital mutual 
information can also be carried out according to the pro- 
cedure described in Sec. II and a block entropy profile 
shown in Fig. [9f b) is obtained. By comparing Figs. IHJb) 
and[5]Jb) it is clear that in the latter case the block en- 
tropy is highly localized, i.e., it takes large values only 
for a few superblock configurations. Repeating the cal- 



culation with fixed M = 64 block states again but using 
an optimized ordering and CAS vector yields a ground 
state energy estimate of -E pcr oxo = —541.55533 Hartree. 
The fast convergence is shown in Fig. [7] by the diamond 
symbol. In fact, the ground state energy after the first 
half sweep, i.e., at the end of the CI-DEAS procedure 
is already far below the one given by Ref. [42| and the 
block entropy profile has the same structure as the one 
obtained after eight sweeps. This can be attributed to 
the CI-DEAS procedure for selecting the appropriate en- 
vironment block states at each iteration step by taking 
care of the renormalized system block and the inclusion 
of the highly entangled orbitals from the beginning. 

There is no need to use additional methods to guaran- 
tee con verg ence like white noise [25j or perturbative cor- 
rection 43j which is an interesting observation in view of 
the results of Kurashige and Yanai. The optimized or- 
dering and CAS vectors used in the calculation are sum- 
marized in Table HI Carrying out the same procedure for 
the bis(/Lt-oxo) isomer of [Cu2 02] 2+ , the convergence of 
the ground state energy with fixed M = 64 states without 
optimized ordering and CAS vector is shown in Fig. llOf a) 
by the square symbol. We also observe for the bis(/it-oxo) 
isomer, that our DMRG energies reproduce the results 
given in Ref. [42f with a fraction of the computational 
resources. It is clear that the convergence of the en- 
ergy is very slow to the value given in Ref. [421 ] . how- 
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TABLE I: Optimized ordering and CAS vector for the [i — rj 2 
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FIG. 11: (Color online) Similar to Fig. [2]but for the bis(/i-oxo) 
isomer of [Cu2 02] 2+ . 



m' < b > 

' < 


> 


! i 

i i 




1 1 

i .. i 


\ f 

ll 


- 

1 sf « , 'i 






11 



50 100 150 200 250 300 

Iteration step 



FIG. 10: (Color online) (a) Convergence of the ground state 
energy in Hartree for various parameters for the bis(/i-oxo) 
isomer. The dashed line gives the ground state energy ob- 
tained in Ref. [i^] using M = 800 block states, (b) Number 
of block states as a function of iteration steps for x = W~ 4 . 



ever, after the entropy-based optimization procedures, a 
much faster convergence to a much lower value can be 
reached. As a result, the optimized CAS vector based 
on the single-orbital entropy shown in Fig. [TT] is given 
in Table ITT] The highly localized two-dimensional entan- 
glement matrix for the optimized ordering is understood 
from the comparison of Figs. Q2] and [T31 The optimized 
ordering vector utilized in the calculation is given in Ta- 
ble In] The entanglement distance idist = 1043 for the 
energetical ordering is reduced to Idist = 102.5 for the 
optimized ordering. The resulting entropy profiles as a 
function of sweeps are shown in Fig. Q31 Again the block 
entropy takes large values only for a few iteration steps 
within a full sweep. 

Since a very important observable related to transi- 
tion metal clusters is the energy difference between the 
two isomers, the ground state energies should be calcu- 



ORD = [ 44 41 43 40 42 31 30 33 32 29 11 10 9 8 7 
20 21 19 18 22 1 2 5 6 4 3 12 17 16 15 
14 13 35 34 36 38 37 39 25 26 24 27 28 23 ] 

CASV = [ 35 13 14 34 3 25 26 18 7 4 8 36 22 15 12 
2 9 10 19 40 5 42 1 11 29 24 43 44 33 41 
32 30 6 21 17 31 27 20 28 38 23 29 16 37 ] 

TABLE II: Optimized ordering and CAS vector for the bis(/i- 
oxo) isomer of [Cu202] 2+ . 
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FIG. 12: (Color online) Similar to Fig. 0] but for the bis(/u- 
oxo) isomer of [Cu202] 2+ . 
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FIG. 13: (Color online) Similar to Fig. [£2 but for the opti- 
mized ordering. 
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lated for the same error margin. This cannot be guar- 
anteed with a fixed number of block states but it can 
be achieved by the DBSS procedure (or an automated 
Richardson- type error protocol [H|). We have thus cal- 
culated the ground state energies of the two isomers us- 
ing the same parameter set, namely for M m j n = 64, 
X = 10~ 4 , Aistart = 256. The convergence of the ground 
state energy is plotted in Fig. [7](a) by the triangle sym- 
bol, while the number of block states selected dynam- 
ically in Fig. [T^b). Since the block entropy is highly 
localized for the optimized ordering for most of the su- 
perblock configurations the number of block states were 
determined by M m i n as can be seen in Fig. [TJb) and 
Fig. ITOTb), We found E peT0 * = -541.55628 Hartree and 
-Ebisoxo = —541.51416 Hartree yielding an energy differ- 
ences of 0.04212 Hartree, i.e., 110 kJ/mol. For both iso- 
mers we have obtained (S 2 ) = 10~ 4 for the ground state 
as expected from the given error margin. 

The convergence of the ground state energy for the 
bis(^i-oxo) isomer is shown in Fig. UOf a) by the trian- 
gle symbol. The fast convergence is again evident. The 
effective size of the environment block during the first 
N iteration steps, i.e., for the CI-DEAS procedure, was 
found to be less than five orbitals. Therefore, using a 
larger number of block states for the environment block 
(M s tart) increases the computational time insignificantly. 
This result indicates that the application of DBSS and 
CI-DEAS procedures in the entropy-based optimization 
guarantees that the DMRG algorithm can be used in a 
black-box fashion, and that chemical accuracy can be 
reached using very limited computational resources. In 
fact, the reduced effective system size in the CI-DEAS 
procedure allows one to obtain the most relevant charac- 
teristics of the entropy functions within a few minutes. 

To ensure that during the CI-DEAS procedure an even 



FIG. 14: (a) Single-orbital, (b) block-entropy, and (c) mutual 
information obtained for bis(/i-oxo) isomer of [Cu202] 2+ as a 
function of DMRG sweeps with optimized ordering and CAS 
vector and for setting the quantum information loss to x — 
10 -4 . The CI-DEAS initialization procedure corresponds to 
symbols with red circle. The optimized ordering and CAS 
vector is given in the text. 



better represented environment blocks are constructed, 
we have repeated the calculations using M m j n = 256 
and M s tart = 1024. After eight sweeps, we obtained the 
following -Eperoxo = —541.57900 Hartree and -Ebisoxo = 
—541.53599 Hartree with a gap of 113 kJ/mol. In addi- 
tion, we have performed more accurate calculations for 
X = 10~ 5 , M min = 256, and M sta rt = 512. The max- 
imum number of block states selected dynamically was 
around 1000 for the peroxo isomer, while for the bis(/i- 
oxo) isomer we found that slightly more states were re- 
quired to reach the same accuracy, so that M max was 
in the range of 1200. In both calculations, we obtained 
(S 2 ) = 10~ 5 as expected for a singlet state. For the given 

= -541.58050 Hartree, 



error margin, we found E, 



peroxo 

-Ebisoxo = —541.53565 Hartree, thus the difference is 
0.04485 Hartree that is 118 kJ/mol. 

In Table IIIIl all relevant DMRG calculations for the 
bis(/x-oxo) and /i — rj 2 : rj 2 peroxo isomers of the binu- 
clear copper cluster are listed. Total electronic energies 
of the two isomers reported by Malmqvist et al. [59| and 
Kurashige and Yanai [43| were calculated in a different 
active space. The latter authors then carried out new 



DMRG calculations 34] , where they reduced their active 
space, but applied the canonical-transformation (CT) ap- 
proach [33| and optimized the orbitals. A comparison of 
absolute energies is therefore not very meaningful as the 
total energy depends on the size of the CAS and the type 
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Ref., method 


^bisoxo -tvpcroxo 


AE 


Reference energies 


[41],CASSCF(16,14) 


-541.50307 -541.50345 


1 


[41J,CASPT2(16,14) 


-542.06208 -542.06435 


6 


[41J,bs-B3LYP 


-544.19419 -544.27844 


221 


[59],RASPT2(24,28) 




120 


Previously published DMRG energies 



[42],DMRG(26,44)[800] -541.46779 -541.49731 78 

[23J,DMRG(26,44)[128] -541.47308 -541.51470 109 

[43J,DMRG(32,62)[2400] -541.96839 -542.02514 149 

[34J,DMRG(28,32)[2048]-SCF -541.76659 -541.80719 107 

[34] ,DMRG(28,32) [2048]-SCF/CT 113 



DMRG energies 


from this work 






DMRG(26,44)[64/256/10~ 4 ] 


-541.51416 - 


-541.55628 


111 


DMRG*(26,44)[256/512/10~ 4 ] 


-541.53499 - 


-541.57896 


115 


DMRG* (26,44) [256/1024/10~ 4 ] 


-541.53599 - 


-541.57900 


113 


DMRG(26,44)[256/512/10 -5 ] 


-541.53565 - 


-541.58050 


118 



TABLE III: The relative energies for the bis(^i-oxo) and 
fi — rj 2 : rf peroxo [Cu202] 2+ isomers obtained in this work 
and from previously published studies are listed in kj/mol. 
All total energies are given in Hartree. The CAS is denoted 
in parentheses as '(electrons,orbitals)' while the information 
on DMRG block states M is given in brackets. Note that total 
energies for different CAS and different types of orbitals can- 
not be directly compared. All results from this work employ 
entropy-based optimization applying the DBSS and CI-DEAS 
procedures. The "*" denotes keeping slightly more DMRG 
environment states during the initialization than necessary 
for a given quantum information loss \- The square brackets 
[M m i n /M a tart/x] state that the DMRG calculation starts with 
M s tart block states, and the minimal number of block states 
is set to M m i n , respectively. 



of orbitals chosen. Instead we shall focus on the relative 
energy which is central to the chemistry of such systems. 
Our entropy-based DMRG calculations agree quantita- 
tively with the RASPT2 (24,28) and DMRG-SCF/CT 
calculations, even though no procedure — such as canon- 
ical transformation — has been employed to account for 
dynamical correlation in our studies. However, this is not 
surprising because Yanai et al. found [34] that the effect 
of CT is only about 4 kJ/mol in the given one-particle 
basis set. In the first study 42] on this problem we al- 
ready found the relative DMRG energies converge faster 
than the absolute energies of each isomer. This finding 
has been confirmed by Yanai et al. [3~i| . 

It is evident from Table IIIII that already our first 
DMRG estimate [42| for this relative energy was much 
closer than the CASSCF result to the correct splitting, 
which we may expect between 110 and 160 kJ/mol based 
on the RASPT2 calculation and the Kurashige- Yanai 
DMRG result for a much larger active space. Our im- 
proved result for the splitting [23j j , which we obtained for 
the original active space but with a reduced number of 
block states Mi = 128 is with 109 kJ/mol already very 



close to the DMRG-SCF/CT and RASPT2 results of 113 
and 120 kJ/mol, respectively. Even the total energies 
turned out to be improved and in the case of the per- 
oxo copper cluster below the small-CAS CASSCF result 
as it should be. In view of the fact that a discrepancy 
in relative energy of about 5 to 10 kJ/mol is acceptable 
for chemical accuracy, we emphasize that both results, 
113 and 118 kJ/mol, are close and that hence already 
the small DMRG calculation is a feasible means to yield 
such relative energies. 



B. Entropic analysis of the two isomers 

Besides the relative energies of the two isomers, one 
might deduce more information from the entropy func- 
tions related to the chemical properties of the binuclear 
copper cores. By comparing Figs. [2] and QTJ it is clear 
that the importance of the orbitals is different for the 
two isomers as shown by the different distribution of the 
single-orbital entropy. In general, almost all orbitals pos- 
sess some 10-20% larger entropy for the bis(/i-oxo) iso- 
mer but some of the orbitals have 2.5-3 times larger 
entropy than those of the peroxo isomer. For example, 
orbitals 3, 14, 34 produce much larger contributions to 
the total entanglement in case of the bis(/x-oxo) isomer. 
The total quantum correlation encoded in the ground 
state, Trot = Ei=i jv s (l)i> i s 3.49 and 5.39 for the 
fx — rj 2 : rj 2 peroxo and bis(/x-oxo) isomers, respectively. 
Thus, the bis(^-oxo) isomer is more entangled which is 
also reflected by the larger maximum value of the block 
entropy for the bis(/^-oxo) isomer as can be seen by com- 
paring Figs. [§] and [T3J 

The highly entangled molecular orbitals for the n — rj 2 : 
rj 2 peroxo and bis(^t-oxo) isomers are shown in Figs. [15] 
and [TjJJ The molecular orbital analysis shows that the 
highly entangled orbitals have an occupation number 
that strongly deviates from either being doubly occupied 
or empty. The entropy-based optimization scheme accu- 
rately determines those orbitals and thus allows one to 
perform efficient DMRG calculations with the smallest 
possible active space for a desired accuracy in the DMRG 
energy. The selected orbitals are close to the Fermi sur- 
face and would be included in a standard CASSCF cal- 
culation if one could employ these large active spaces as 
in DMRG calculations. 



V. SUMMARY AND PERSPECTIVE 
A. Relation to Tensor Network States 

Our procedure has been demonstrated on a one- 
dimensional spatial topology related to the DMRG 
method, while the obtained two-dimensional entangle- 
ment matrix could be employed more efficiently in meth- 
ods based on tensor network state (TNS) approaches. For 
example, the two-orbital mutual information can provide 
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FIG. 15: (Color online) Molecular orbital pictures of the highly entangled orbitals for the bis(/i-oxo) isomer. The orbitals were 
selected with respect to their one-site entropy as shown in Fig. 1111 The number below each orbital corresponds to the orbital 
index and the occupation number is written in the parentheses. 
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FIG. 16: (Color online) Molecular orbital pictures of the highly entangled orbitals for the peroxo isomer. The orbitals were 
selected with respect to their one-site entropy as shown in Fig. [2] The number below each orbital corresponds to the orbital 
index and the occupation number is written in the parentheses. 



an optimal value for the coordination number of each or- 
bital in case of the tree tensor network state (TTNS) 
approach [Tll j and an optimized network topology can 
be determined. Since in this latter method the distance 
between highly entangled orbitals can be reduced signif- 
icantly compared to the DMRG topology entropy-based 
optimizations are expected to boost its convergence prop- 
erties to a great extent. As an example, we present cal- 
culations performed with the TTNS method on the Be 
atom studied recently in Ref. [ll|. The components of 
the two-orbital mutual information larger than 10 -4 are 
shown in Fig. ITTTa) . The convergence of the ground state 



energy with a fixed bond dimension D = 2 and coordi- 
nation number z — 3 for three different topologies shown 
in Fig. HU is plotted in Fig. Tf\h). It is found that 
for this set-up the relative error in the ground state en- 
ergy dropped by an order of magnitude if an optimized 
topology based on the two-orbital mutual information is 
used. For the one-dimensional DMRG topology and for 
the energetical ordering, we found that Idist = 0.7251 
whereas it is reduced to 0.236 for the optimized topol- 
ogy of the TTNS method. For the topologies shown in 
Fig. [TS] (a) and (b) we found a considerably larger val- 
ues. We want to emphasize that for the one-dimensional 
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FIG. 17: (Color online) (a) Similar to Fig. Q] but for the Be 
atom, (b) Convergence of the ground state energy for the 
Be atom obtained with the TTNS method with a fixed bond 
dimension D — 2 and coordination number z = 3 for three 
different network topologies shown in Fig. [18] 60] . 
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FIG. 18: (Color online) Three different network topologies 
were used for the TTNS method. The topology shown in (c) 
is optimized based on the two-orbital mutual information. 



DMRG topology the lowest value of Tdist is 0.30 by in- 
cluding all permutations of the orbitals. For the opti- 
mized ordering based on the Fiedler vector, Idist is 0.319 
for which, in fact, we obtained the lowest ground state 
energy (E rel = 1.69 x 10~ 5 ) when fixed M = 8 DMRG 
block states were used. A detailed study based on spec- 
tral analysis and graph theory and extension of our ap- 
proach to larger systems is under progress and will be 
part of a subsequent work. 

The optimization of the spatial arrangement of or- 
bitals and network topologies has a significant influence 
on the convergence properties of MPS and TNS algo- 
rithms. However, the total entanglement encoded in the 
wave function, Trot, cannot be changed. The optimiza- 
tion of the one-particle basis can yield entanglement re- 
duction in the system and the value of Trot can be ma- 
nipulated, consequently. In fact, as discussed in Ref. [l5[ 
by constructing an optimal basis it might be that order- 
ing is either obvious or irrelevant. In the past few years, 



various procedures have been developed for orbital op- 
timizations [HI HH, I3214371 1 . We will therefore continue 
to examine mutual information with respect to orbital 
optimization in our subsequent work. 



B. Conclusion 

In this paper, we have studied a transition metal clus- 
ter from a quantum information theory perspective using 
the DMRG method. By calculating various entropy func- 
tions and the two-orbital mutual information we have 
proposed recipes to perform DMRG calculations in a 
black-box fashion on complex chemical compounds. 

Optimizations based on the two-orbital mutual infor- 
mation can be related to graph theory and spectral anal- 
ysis of seriation problems. Our cost function is inter- 
preted in terms of the Fiedler vector of the graph Lapla- 
cian. Our results confirm the importance of taking en- 
tanglement among molecular orbitals into account and 
the usefulness of graph theory for carrying out efficient 
calculations. Chemical characteristics of the two isomers 
of [Cu202] 2+ have also been analyzed and interpreted in 
terms of the entropy functions. 

The present work confirms our previous findings that 
even small- M DMRG calculations provide a qualitatively 
correct description of transition metal clusters as demon- 
strated in Ref. [Hill]. 

Our results indicate that highly entangled orbitals 
form subgroups. Therefore, a coarse-graining approach 
might be possible which could be efficiently implemented 
by the multiscale-entanglement-renormalization ansatz 
(MERA) J|. 

In future work, we shall explore the promising 
DBSS/CI-DEAS method in DMRG calculations on a set 
of test molecules featuring states whose relative energy 
is difficult to calculate. Also, additional options for im- 
provement like extrapolation schemes [lij shall be inves- 
tigated. 
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